Density pattern in supercritical flow of liquid ^He 
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A density functional theory is used to investigate the instability arising in superfluid *He as it flows 
at velocity it just above the Landau critical velocity of rotons Vc- Confirming an early theoretical 
prediction by one of us [JETP Lett. 39, 511 (1984)], we find that a stationary periodic modulation 
of the density occurs, with amplitude proportional to (u — Vc)^^^ and wave vector equal to the roton 
wave vector. This density pattern is studied for supercritical flow both in bulk helium and in a 
channel of nanometer cross-section. 

PACS numbers: 68.10.-m , 68.45.-v , 68.45.Gd 



According to Landau criterion superfluid motion of 
^He moving with velocity u is possible only if m < v^, 
where the critical velocity Vc is equal to the slope of the 
tangent to the roton part of the spectrum {vc ~ 56 m/s). 
As the Landau critical velocity Vc is reached, the liquid 
becomes unstable against a spontaneous excitation of ro- 
tons. Reaching the roton critical velocity is difficult in 
practice since other type of excitations, e.g. quantized 
vortices, are produced in bulk "^He well below Vc- How- 
ever, the occurrence of vorticity could be suppressed by 
allowing ^He to flow through very narrow channels. In 
fact, the critical velocity for the creation of vortex pairs 
in a channel of diameter D is t;™'-*^^ {h/MD) ln(Z)/^) 
0, where ^ is the '''He healing length ~ 1 A), so that 
it can exceed the roton critical velocity for channels of 
nanometer size J^. 

Several years ago a theoretical prediction was made by 
one of us ^3] that superfluid ^He flowing at a velocity 
greater than the Landau critical velocity Vc should un- 
dergo a phase transition from a spatially homogeneous 
state to a layered state characterized by a periodic den- 
sity modulation in the direction of motion. Such a modu- 
lation is stationary in the frame moving with the fluid and 
has a characteristic wavelength A = 2ttH/pc ~ 3.58 A, 
where Pc is the roton momentum. This prediction was 
derived within a simplifled model describing a weakly in- 
teracting roton gas with coupling constant g. The nature 
of the transition was found to depend on the sign of g: if 
g > {g < 0) the transition is predicted to be continuous 
(discontinuous). In Ref. the estimate g = 2 x 10""^^ 
erg cm'^ jj| was used and the amplitude of the density 
modulations was found to be @ 
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where po is the bulk density and |ylp(5(/ia; — e{pc)) is 
the roton contribution to the dynamic structure factor 
S(q,uj). In Ref. the latter was estimated by ignor- 
ing the multiphonon part of S{q, uj) and using the /-sum 
rule. A better estimate can be extracted from neutron 



scattering experiments p, where one finds |Ap ~ 0.9. 
Inserting this value in Eq. one gets 



Ap/po^3[{u-v,)/v,]'/^ . 



(2) 



The occurrence of this stationary nonuniform state 
originates from the presence of a pronounced minimum 
aXp — Pc in the bulk "^He spectrum, e(p), and is similar to 
the structural phase transition in crystals induced by the 
softening of phonon frequencies with some defined wave- 
length. It is also worth mentioning that similar periodic 
modulations in the ^He density profile near a moving im- 
purity have been observed in recent computer simulations 
[3, although not explained in terms of this instability. 

In this work we investigate the occurrence of density 
patterns in the supercritical "^He flow by performing den- 
sity functional (DF) calculations. We consider a uni- 
form flow in bulk liquid (with no vorticity) as well as 
in a nanochannel. We use the DF approach proposed 
in Ref. 8] and later improved in Ref. 9], which gives a 
quite accurate description of inhomogeneous configura- 
tions of liquid *He at T = 0. The energy of the system 
is expressed as: 



EM ^J^l dri^Vpf + I dvEc{v) 



(3) 



The explicit form of the correlation energy is given in 
Ref. ISj. The static equilibrium profile p(r) in an arbi- 
trary external potential can be obtained by minimizing 
the functional Eq[p] with respect to density variations, 
subject to the constraint of a constant number of atoms. 
The dynamics can be studied as well by means of the time 
dependent DF method, with the DF proposed in Ref. 
playing the role of the effective Hamiltonian driving the 
time evolution of the system. In the dynamical case, 
the functional contains an explicit dependence upon the 
local current density field j (r) through a phenomenologi- 
cal term which accounts not only for the usual hydrody- 
namic current density but also for non-local "backflow" 
effects. The resulting DF (named Orsay- Trento Func- 
tional), which will be used in our calculations, has the 
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following form: 

E[p,^,]=EQ[p\+ J dvHj . (4) 

An appealing feature of the above functional, which turns 
out to be essential to perform accurate time dependent 
DF calculations is that it reproduces quantita- 

tively not only a number of static properties, but also 
the observed phonon-roton spectrum of bulk "^He. 

The minimization of the above density-current func- 
tional, subject to the constraint of a fixed number of 
*He atoms and of fixed total momentum, can be done 
in practice by evolving in the imaginary time domain 
a non-linear Schrodinger equation for the order param- 
eter vE'(r), where the Hamiltonian operator is given by 
H = -1/(2M)V2 + C/[p,v]. The effective potential U is 
defined in terms of the variational derivative of the en- 
ergy functional, and its explicit expression can be found 
in Ref. jlfl] . From the knowledge of the complex wave- 
function ^ = (/)e''^ one can get immediately the density 
p(r) = (j)^ and the velocity field v(r) = (1/M)V6. Since 
we are interested in stationary states of ^He in the pres- 
ence of a uniform flow, we minimize the above functional 
in the frame of reference moving with the liquid, which 
we assume to flow with some given velocity u along the 
X-axis: The Hamiltonian density H thus acquires an ad- 
ditional term H' — H — uPx, Px being the ''He total 
momentum component along the direction of motion. 
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FIG. 1: Amplitude of the density modulation along the di- 
rection of *He motion, computed during the functional min- 
imization. Solid line: u — 1.29 Vc, dotted line: u = 1.14 Uc, 
dashed line: u = 0.99 Uc, dash-dot line: u — 0.84 Wc 

First, we address the problem of the Landau roton in- 
stability in bulk. As discussed above, we expect that 
when u > Vc the uniform density configuration is not sta- 
ble, but it is instead a metastable state corresponding to a 
saddle point of the energy landscape of ''He. In our case, 
the systems is allowed to reach the lowest energy con- 
figuration by following the (dissipative) imaginary-time 
evolution. The calculation is performed in a periodically 
repeated supercell where the size of the cell along the x- 
direction (which we take as the direction of "^He motion) 



is L. Our procedure to trigger the instability is the fol- 
lowing: we start with the uniform system in the moving 
frame of reference and slightly perturb the (uniform) den- 
sity with a sinusoidal modulation with a small arbitrary 
amplitude and with a wavelength A allowed by the pe- 
riodic boundary conditions in L. We then minimize the 
functional in the frame of reference moving with some 
chosen velocity u, with the only constraint of a constant 
number of ''He atoms. If i or A are not a multiple of the 
characteristic wavelength Ac = (2'K)/kc {kc = Pc/fi being 
the Landau critical wave- vector) then, irrespective of the 
initial perturbation and of the particular value chosen 
for u, the perturbing modulation rapidly smoothes out 
during the minimization, and the uniform liquid state is 
recovered as the minimum energy configuration. 
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FIG. 2: Density profiles along the direction of ^He flow x-axis. 
The three profiles have been calculated, in order of increasing 
amplitude, with u = 1.07 , 1.22 , 1.37 «c, respectively. 

The Landau instability shows up when the size of the 
cell is such to accomodate an integer number of charac- 
teristic wavelengths {L = 7Ac). In this case we indeed 
find that there exists a threshold velocity Vc separating 
two regimes. If w > Vc, the stationary state is charac- 
terized by a density with a periodic modulation of wave- 
length Ac and with an amplitude depending on u. On the 
contrary, when u < Vc the initial modulation is rapidly 
smeared out during the minimization, and again one finds 
that the density of the stationary state is uniform. The 
critical velocity is found to be Vc ~ 58 m/s, which co- 
incides with the minimum value of e{p)/p predicted by 
the same DF and is also very close to the value of the 
Landau critical velocity of rotons as obtained from the 
experimental phonon-roton spectrum. This behavior is 
summarized in Fig. 1, where we plot the evolution of the 
amplitude of the density modulation as it varies during 
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the minimization procedure, for four different ^He veloc- 
ities: the two upper hnes have u > Vc, whereas the two 
lower lines have u < Vc- Note the critical slowing down 
for values of the ''He velocity close to the critical value Vc, 
where a very long imaginary-time evolution is required to 
converge towards the equilibrium stationary state. 

Different stationary density profiles along the direc- 
tion of ^He motion, corresponding to different values 
of u > Vc, are shown in Fig. 2. The average value 
of each curve corresponds to the saturation density of 
bulk ^He, po = 0.0218 A^^. A fit to the calcu- 
lated points shows that their shapes, at least for val- 
ues of u not too large, is almost exactly sinusoidal, i.e, 
p{x) — pq[1 + {Ap/pQ)sm{kcx)]. In Fig. 3 we also 
show the x-component of the calculated ^He velocity 
v(r) = (1/M)V0, in units of Vc, for the same states of 
Fig. 2. Note the oscillating character of the velocity, in 
phase with the density modulation, and the large ampli- 
tude of oscillations, which becomes more asymmetric as 
the velocity increases. The spatial average of the velocity 
profiles shown in Fig. 3 is zero, as expected. 



2 - 




X (A) 



0.8 




FIG. 4: Amplitude of the density modulation a function of 
the fluid velocity. Points: DF results. Line: fitting function 
lM[{u-v,)/v,]^^^ 



Now we investigate the motion of liquid ^He in a nar- 
row channel of nanometer transverse dimensions. We 
consider liquid "^He confined between two infinitely ex- 
tended, weakly attractive planar surfaces separated by a 
very small distance, ~ 50 A. We model the two surfaces 
with an external potential which mimics the adsorption 
properties of the Rb surface, which is the weakest surface 
which is wet at T = by liquid ^He llj. The number 
of ^He atoms in the system is chosen in such a way that, 
when the "^He is at rest, the equilibrium density near the 
center of the channel reaches the value corresponding to 
the saturation density of bulk ^He, po. 



FIG. 3: Velocity profile along the direction of ''He flow a:-axis. 
Same values of u as in Fig. 2 

The main result of this work is summarized in Fig. 4 
where we show the behavior of the amplitude Ap/po for 
u > Vc- We find that the law Ap/po = l-01[(w— Wc)/wc]^^^ 
(solid line) very nicely fits the numerical results (points). 
The velocity dependence is thus the same as in Eq.Q 
except for the different numerical coefficient. Our DF 
calculations are consistent with a repulsive (positive g) 
roton-roton interaction. Using Eq. and the fitting 
coefficient 1. 01, we find g ~ 1.8 x 10"'^^ erg cm'^. It is 
worth stressing that direct measurements of g are not 
available and previous theoretical estimates significa ntly 
differ both in magnitude and sign (see, for instance, 
and references therein). 




FIG. 5: Density profile across the channel section 
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In Fig. 5 we show the density profile along the z- 
direction, i.e. across the channel, for u = 0. The ''^He 
density decreases rapidly to zero near the solid surfaces 
on both sides of the channel due to the ^He-Rb inter- 
action. The same interaction is also responsible for the 
density oscillations near the walls. The dotted line shows 
the value of the bulk saturation density po- Fig- 6 shows 
a contour plot of the density in the xz-plane for the sta- 
tionary state developed at u = 1.22 Vc- The complex 
pattern near the walls is again due to the "^He-Rb in- 
teraction. However, the dominant feature is the density 
modulation along x in the central part of the channel. 
This sinusoidal oscillation coincides, for the same value 
of u, with the one that we already obtained in bulk "^He. 
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ical flow just above Vc and in the absence of vorticity. 
Due to the short wavelength of the density modulations, 
of the order of the atomic spacing, its direct observa- 
tion, with X-rays for instance, might be difficult. Indi- 
rect evidences of the density modulations could however 
be measurable, for example throught their possible ef- 
fects on transport properties. Recently, He adsorption 
within a reg ular porous medium called FSM-16, has been 
studied This silica-based material is characterized 

by ordered arrays of long, uniform pores, with diameters 
ranging from 1.5 to 10 nm. When ^Hc is adsorbed within 
the pores, 1-2 solid-like layers are expected to form, coat- 
ing the internal walls of the pores, leaving however room 
for additional ^He in the liquid state. A pressure gra- 
dient between two open ends of an array of pores could 
in principle be used to force liquid "^He to move through 
this system, until it is expelled from the pore end. If dur- 
ing this process the critical velocity is reached, then the 
occurrence of the above described density pattern might 
induce the fragmentation of the ejected liquid filament 
into regularly distributed nanodroplets. A similar pro- 
cess might occur in the experiments of Ref. j^, where 
liquid ^He is discharged into vacuum through a microm- 
eter nozzle. The structure of the ejected filament was 
interpreted in terms of a Rayleigh instability, but the oc- 
currence of density patterns near the nozzle could also 
play a role 16]. In this perspective, the effects of den- 
sity modulations in "^He supercritical flow in this type 
of experiments deserve further investigations. Finally, it 
is worth mentioning that a similar phenomenon may oc- 
cur in Bose-Einstein condensed gases with dipole-dipole 
interactions [l^ . 



FIG. 6: Contour plots of the density along the channel 

Our DF calculations support the predictions of Ref. Q 
on the occurrence of a density pattern in the supercrit- 
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